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Abstract. An efficient microcanonical dynamics has been recently introduced 
for Ising spin models embedded in a generic connected graph even in the presence 
of disorder i.e. with the spin couplings chosen from a random distribution. Such 
a dynamics allows a coherent definition of local temperatures also when open 
boundaries are coupled to thermostats, imposing an energy flow. Within this 
framework, here we introduce a consistent deflnition for local energy currents and 
we study their dependence on the disorder. In the linear response regime, when the 
global gradient between thermostats is small, we also define local conductivities 
following a Fourier dicretized picture. Then, we work out a linearized "mean- 
fleld approximation" , where local conductivities are supposed to depend on local 
couplings and temperatures only. We compare the approximated currents with 
the exact results of the nonlinear system, showing the reliability range of the 
mean-fleld approach, which proves very good at high temperatures and not so 
efficient in the critical region. In the numerical studies we focus on the disordered 
cylinder but our results could be extended to an arbitrary, disordered spin model 
on a generic discrete structures. 

keywords: Transport processes (Theory), Heat conduction. Disordered systems 
(Theory) 



1. Introduction 

The transport properties featured by systems in stationary states, far from equilibrium 
is of both theoretical and practical interest: On the one hand, there exist non trivial 
problems (e.g. validity bounds of the Fourier description at microscopic scale, influence 
of spatial and topological inhomogeneities featured by the substrate, etc.), which still 
lack an exhaustive solution; on the other hand, from nanoscales to biological matter, 
the emergence of new materials poses a challenging number of practical problems 
(e.g. fluctuations, clusterizations, correlations, role of geometrical irregularities, etc.), 
where the previous, theoretical approaches get an applicative relevance [U [H |3l |4]. 
Discrete models, from simple or interactive random walks to classical spin models, 
play an important role for all such questions [5, 6j ^ HI |9l [ini HU [I^- In particular, 
here we deal with an Ising system coupled with thermostats imposing an energy flow, 
and we study its behavior at microscopic scales. 
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The dynamics we adopted in order to simulate the system evolution is a 
microcanonical dynamics, recently introduced in [5], which has allowed for important 
enhancements in the study of heat flow in a spin lattice coupled to thermostats at 
arbitrary fixed temperatures. The advantages of this dynamics include the simple 
definition of a local temperature on each link, the capability of a direct computation 
of the conductivity (instead of mere difFusivity) , the possibility of applying to an 
ample variety of discrete spin systems ranging from models with random couplings 
to models defined on inhomogeneous networks [13]. In particular, the cases of 2- 
dimensional cylindrical systems, with ordered or random distributions of couplings 
Jij > 0, have been considered in [6]. The cylindrical geometry of the lattice suggested 
to average observables (temperature, magnetic energy, etc.) not only in time but also 
along the "columns" , i.e. periodic rings at equal distance from the thermostats, and 
hence orthogonal to the flux flowing along "rows" . This procedure clearly accelerates 
the numerical convergence to well stabilized stationary values. In this frame, the 
validity of the Fourier picture for the energy transport is ensured only on the average, 
by verifying that the average energy flux through a column is proportional to the 
average temperature difference between the sides of the column. However it is clear 
that, due to quenched disorder, the local temperature fluctuates even within a single 
column. A natural question is then whether such spatial fluctuations influence the 
local transport properties. 

In this perspective, a basic problem, constituting a non trivial task in itself, is 
a self-consistent definition of microscopic currents able to account for the substrate 
inhomogeneity (here generically referring to topological disorder and/or non-constant 
couplings). In the present paper we give a consistent solution to this problem which, 
remarkably, applies to arbitrary, connected supports. The microscopic link currents 
originating in each microscopic move can be defined indeed independently of the 
subjacent geometry. 

Afterwards, resuming the case of the disordered ferromagnetic cylinder, this very 
general definition is used to check the validity of the Fourier law at microscopic scale. 
First, we verify to this end that, locally, the currents depend on the disorder in a 
non trivial way. More precisely, local currents and temperatures depend not only on 
the relative coupling Jij, but also on the whole configuration of magnetic couplings. 
Second, we show that, notwithstanding such a complex behavior, a local conductivity 
can be defined on each link. In fact, in a linear response regime (i.e. for small enough 
temperature differences AT at the borders), the ratio between the local currents and 
the local temperature gradients are independent of AT, providing a good definition 
for the local conductivities Kij. Clearly, each Kij also depends in a non trivial way 
on the system parameters, its value being indeed determined not only by the local 
coupling and temperature but also by the actual realization of the global disorder. 

Despite such a complex situation, a simple correlation between local conductivities 
and local couplings is present, since on the average the conductivities increase with 
Jij. Hence, we introduced an approximated linear fitting Kfj^ = AJij + B, where 
Kfj^ are the approximated "mean field" conductivities and the constants A and B 
depend on the system parameters only, i.e. they are independent of the disorder 
realization. Within this approach the conduction properties (i.e. the currents) can 
be simply evaluated by solving the discretized linear Fourier equations with suitable 
border conditions. The approximated currents obtained in this way may be compared 
with the exact results of the fully non linear system. We evidence that the approach 
works quite well at large enough temperatures, while near Tc (the critical temperature 
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of the equilibrium Ising lattice) correlations between local conductivities and currents 
play a much more important role, and the local "mean field" approach becomes less 
efficient. This linear approach for the evaluation of local currents allows for example 
to detect links which, due to the particular disorder realization, are expected to be 
traversed by large currents and therefore to locate the possible critical bonds where 
the system may fail because of too large a load. 

The paper is organized as follows: In Sec. [2] we briefly review the microcanonical 
dynamics exploited to simulate the evolution of the system and in Sec. [3] we introduce 
a consistent definition of local currents. Then, in Sec. S] we present the numerical 
results focusing on local currents, on local conductivities and on their dependence on 
the coupling pattern, while in Sec.[5]we show a mean-field approach to get an estimate 
of local currents as a function of a given, arbitrary set of couplings; finally Sec. [S] is 
left for conclusions and perspectives. 

2. Microcanonical Dynamics 

The very novelty of the dynamics introduced in [5] consists in assigning to the links, 
beside the usual magnetic energy i?™ due to the spin configuration, a "kinetic" energy 
E^j , i.e. a non- negative definite quantity whose variation can compensate the positive 
or negative gaps of magnetic energy determined by the spin flips of the adjacent nodes. 
This is clearly different from the Creutz microcanonical procedure [Mj [15] , where the 
required energy compensation is extracted or assigned from a bounded amount of 
energy lying on the nodes themselves: in the latter case, in fact, both the energy 
boundness and the connectivity of the nodes could entail many limitations, ranging 
from geometrical or topological constraints to the non ergodicity of the system at low 
energy density. 

Precisely, an elementary move consists in the following: 

(i) starting from a random distribution of link energies, extract randomly a link (i, j) 
and one of the possible four spin configurations for it; 

(ii) evaluate the variation AE'™ of the magnetic energy due to this choice, checking 
the whole neighborhood of the link. Obviously, if couplings are ferromagnetic, 
being Jnk € [1 — e , 1 + e] for any link (n, k) and e > 0, then A£'™ is a real number, 
and if e = then Ai?'" is an integer whose value depends on the connectivity of 
the k and n nodes; 

(iii) if AE™ < 0, accept the choice and increase the link kinetic energy E^j of AE"^; 

(iv) if AE"^ > 0, accept the choice and decrease the link kinetic energy of Ai?'" only 
if the link energy remains non negative. 

The unit time step will be a series of moves, where N is the number of links. 
Since a link is defined only by the adjacent nodes, and nothing in the rule above refers 
to a definite structure (e.g. a lattice), the move is implementable on every non-oriented 
connected graph. This ensures the great generality of the dynamics. 

Starting from the observation that magnetic and kinetic energies behave as 
non correlated observables, in [6J many points have been supported by theoretical 
arguments and tested numerically, for both homogeneous and disordered links. In 
particular, at equilibrium the Boltzmann distribution is recovered and the system is 
ergodic at all temperatures. This opens the possibility to extend the definition of 
temperature Tij as a local observable, in fact a link observable, which recovers the 
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averaged kinetic energy (E^j). This holds also in non-equilibrium, stationary states, 
i.e. states forced by thermostats at different temperature. 

As for thermostats, their definition requires in general the presence of contact 
borders, which for the cylinder are the first and last columns. In this case, every 
thermostat consists in a number of additional columns (2 are enough) regulated at 
every step by the usual Metropolis equilibrium dynamics with the wanted inverse 
temperature /3. 



3. Definition of Local Currents 



In order to calculate local conductivities also in the presence of disorder (due to either 
topological and coupling inhomogeneity) it is necessary to define a current for each 
link. We now introduce the scheme through which we are able to consistently assign 
a current to any arbitrary link in the structure considered. Such a scheme can be 
applied to a generic structure, as envisaged in Fig. [T] 

In our dynamics energies are naturally assigned to each links; however in order 
to define link currents it is useful to assign one half of the link energy to each of its 
relevant nodes. Let us consider the link (i, j) connecting i and j; Let Vi and Vj their 
neighborhoods, Zi =\ Vi \ and zj =1 T^" | their respective coordinations. Then, k G Vi 
and h eVj are the neighbouring site labels for i and j, respectively. We fix a direction 
for currents on adjacent links: currents on links connected through i are incoming 
while those connected through j are outgoing. 

Now, given a spin flip involving the link under consideration, namely either the 
i-th spin, or the j-th spin or both, for adjacent links as well as for {i,j) we have a 
possible energy variation denoted with AiJ^fc, AEjh and AEij respectively. For links 
(«, k) and {h,j) half of such variations contributes to the current on the link itself since 
it represents the energy flow outgoing and incoming from the sites k and h. Therefore 
we have 

J _ AEki _ AEjh 

J-ki — ^ , Jj/i — H ^ , 

where different signs derives from the flow direction we have chosen (see Fig. [1]). 
As for lij and the energy variation AEij , they satisfy the following 

fceVi k£Vi 

In the above equations, the left-hand side represents the currents arriving and 
departing from i (and j) while the right-hand side is the consequent energy variation. 
Hence, we get AEij = ^'llkev- ^-^ki — J2hev-^-^hj^ consistent with energy 
conservation, and 



- E ^Ek^ + E ^Ejh 



(3) 



This scheme works for any arbitrary topology and, of course, even in the presence 
of a disordered distribution of couplings: It only requires the knowledge of the local 
energy variations consequent to any spin-flip. 
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Figure 1. Schematic representation of currents for an arbitrary link belonging 
to an arbitrary graph. In this example the link considered is the one labelled as 
with Zi = A while Zj = 3. The arrows indicates the direction attributed to 
the energy flow on each link. 



4. Numerical results 

In the following we report and discuss the results obtained by means of Monte Carlo 
simulations performed on squared cylinders oi N = LxL sites, for different realizations 
of disorder (encoded by the N x N matrix J). The first and last columns of L sites are 
open, in contact with thermostats at temperatures Ti and T2 (to fix ideas Ti < T2), 
and different choices of temperatures are considered. Since we are interested in local 
quantities, we especially focus on small sizes, which allow fast thermalization though 
displaying the relevant features of the non-equilibrium behavior (see [5]). 

First of all, let us consider local currents and local temperatures. For a single 
realization of the disorder in the window (1 — e, 1 + e), currents lij may be calculated 
according to the scheme described in the previous section; one can also measure the 
local temperatures Ty which, for randomly distributed couplings equal the relevant 
average kinetic energy (E^^) [B]. Then, from such local temperatures Ty, we can 
estimate a temperature Ti to associate to each node, namely 



Therefore, the local temperature gradient among nodes i and j is naturally given by 



Numerical data for the average currents (lij) and ATij, as a function of the 
pertaining coupling strength Jij, are shown in Fig. [2] similar results are obtained for 
different realizations J. We notice that, being AT = T2 — Ti the global difference of 
temperature, local gradients are distributed around the expected value AT/L, with 
large spread especially for low temperatures Ti , T2 and a slight correlation with the 




(4) 
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Figure 2. Local currents lij (left panel) and local temperature gradients AT^j 
(right panel) versus the corresponding coupling Jij for a cylindrical lattice with 
linear size L = 10, AT = 0.2 and different choices for T2, as shown by the legend. 
The coupling pattern J is the same for all sets of data points depicted. 



pertaining couplings, that is, large interaction strengths Jij correspond to smaller 
gradients ATy . Local currents display a larger degree of correlation: large interaction 
strength Jij correspond to larger magnitudes for currents {lij). We also notice that 
different temperatures for thermostats give rise to similar, though shifted, distributions 
of data points. Other realizations of the same disorder give, of course, different point 
distributions, but the linear interpolation and the value of the fluctuations prove to be 
very robust, so that the definition of interpolating currents lij is reliable at every fixed 
e and Ti, T2. ft is also noteworthy that, for a fixed AT, currents are not monotonic in 
T: referring to Fig. [2] (left panel), their magnitude is maximum at T = 3, i.e. around 
the critical temperature expected for the (disordered) two-dimensional Ising model 

m- 

Moreover, when AT ^ I, the linear response theory holds: Both the local currents 
(/y ) and the local temperature gradients ATij are proportional to the global difference 
of temperature AT, as corroborated by the collapse of data points in Fig.[3l where the 
values relevant to different temperatures are compatible with the numerical error. In 
this perspective, given local currents and local gradients, we can introduce the local 
conductivities according to 



Krj = TfT^fT, (5) 
J J J-j 

which are well defined quantities describing the microscopic conduction of the system. 
We expect that a similar definition works as well for generic topologies, at least in the 
regime of small gradients. 

From energy conservation and equation ([5]) we obtain the local Fourier equation 
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Figure 3. Local currents lij divided by the global temperature difference AT 
as a function of the pertaining coupling Jij . As explained by the legend different 
gradients, namely AT = 0.05, AT = 0.1 and AT = 0.2 are compared, while 
T2 = 5, the coupling pattern J and L = 10 are kept fixed. Notice that the 
average difference between data pertaining to AT = 0.05 and AT = 0.1, and data 
pertaining to AT = 0.1, AT = 0.2 is, respectively, 0.3 and 0.4, therefore below 
the error ^ 4%. 



characterized by link dependent conductivities 

j:k.at.-t,)^-^, (6) 

where (Ei) = i/^J^jev i^ij) is the total energy relevant to site i; the continuum 
notation is used for convenience, with the usual warning about the meaning of 
derivatives in these discrete-time systems (see for instance 16 ). In general terms, 
the expression in Eq. |5] describes a system where an external field, or gradient, along 
one axis and a fluctuating local field, or disorder, have been applied; the former makes 
the temperature increase by a constant amount per row of nodes, while the latter 
gives rise to currents non-trivially depending on the whole environment. Indeed, the 
same equation is also used in the context of random resistor networks T7], where the 
voltage and the conductance play the role of the temperature and of the conductivity, 
respectively. 

We remark that conductivities Kij depend on system parameters in a very 
complex way, indeed their values is determined by the temperature T, by the degree 
of disorder e and by the whole coupling pattern J. Moreover, we verified that the 
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Figure 4. Local conductivities Kij versus pertaining local couplings Jij for a 
square cylinder with linear size L = 10 and thermostats at temperatures Ti = 4.8 
and T2 = 5.0, respectively; different degrees of disorder e = 0.05, e = 0.1 and 
€ = 0.2, have been considered and represented with different symbols, as shown 
in legend. Linear curves represent the best fits, whose angular coefficients are 
approximately 358, 279 and 226 respectively. 



dependences on the three arguments are intrinsically interplaying, namely that given 
Kij — /(J, T, e), factorizations like /(J, T, e) — /i(J) • f2{T, e) are ruled out. 

As shown in Fig. 21 local conductivities are correlated on the average with the 
local couplings, and such a correlation gets stronger (the fitting curve has larger slope) 
for smaller values of e. The "ordered system limit" , i.e. e — > 0, is in a sense singular, 
since the distributions of the Jij shrinks in a single point. 

Now, from such local conductivities it is possible to derive an estimate for the 
conductivity K{T, J) expected for a system at a temperature T and in the presence 
of disorder J, by assuming AT <C 1 and averaging over all local conductivities: 

In Fig. [5] we compare such measures realized at different temperatures with a 
"mesoscopic" measure of K based on the heat flow passing from one layer to the 
next one in a similar cylinder . The very good agreement between the two estimates 
provides a further confirmation about the consistency of our definition of local currents 
and conductivities. 
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Figure 5. Comparison between average local conductivities measured according 
to Eq. [7](») and according to a "mesoscopic" measure of conductivity (continuous 
line) 0, respectively, as a function of temperature. 



5. Mean Field Approach 

Within the hnear response approach, local temperatures can be evaluated by solving 
the Fourier equation ^ imposing the stationarity of local energy at every node: 

J2k,j{T,-T,)^0. (8) 

Suitable boundary conditions should be chosen forcing the temperatures at the borders 
to be fixed at T and T + AT respectively. Then currents can be evaluated from 
temperatures as (lij) = Kij{Ti — Tj). Clearly, Eq.[8]is useless for practical calculations, 
as Kij has to be obtained from the (numerical) solution of the whole spin dynamics. 
Moreover, we have already evidenced that Kij^s depend on the whole configuration of 
local couplings Jij in a non trivial way. However Fig. |4] suggests that a simple "mean 
field" approximation should be possible imposing Kij to be dependent only on the local 
coupling Jij, i.e. the mean-field, local conductivities are defined as Kfj^ — AJij + B, 
where A and B are the fitting parameters used in Fig. |4] depending on T and e 
only. Clearly, this approach is much simpler since, once A and B are known, the 
mean-field, local conductivities can be inferred for any realization of the disorder. 
Then, using the linear equations (|5]), one can obtain mean- field, local temperatures 
E,GV, K^^iT^^ - 7f ^) = and currents /^/^ = K^^{T,'^'^ - Tf ^). 
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In Figure [5] we evidence that local, mean-field currents are indeed a good 

approximation of the exact results. More precisely, once the local deviations and 
are defined as 

{I^J)-I^J, (9) 

.MF_jMF_j^^^ (10) 



/MF _ 7-MF f 
fij — J-ij -'iji 

we can quantify the correlation between the real values (lij) of currents, i.e. those 



7-MF 
Hj J 

from the mean-filed approach, by means of the correlation coefficient 



obtained from numerical simulations, and the estimate values Iij , i.e. those obtained 



^_ E(0,,-.0^^)-E,^.,-E0^^ ^^^^ 



y'[E(0Mf)2_(IE^MF)2] 



where averages E are obtained summing over the whole set of links and lij is the 
value of the local current obtained with the linear fit of Figure [2j Notice that C 
ranges from —1 (anticorrelation) to +1 (correlation) and C = means no correlation. 
We measured the quantity C finding strictly positive values for all the temperatures 
considered. The positivity of correlation evidences that if^^ approximates currents 
{lij ) better than lij . This remarkable property is easily explained since the mean field 
field approach not only takes into accounts the correlations between local couplings, 
currents and conductivities evidenced in Figs.j^Hl but also takes into account the local 
conservation of energies encoded in Eq. |8] and representing one of the basic features 
of the microscopic spin dynamics. In particular, on the contrary of lij, currents ifj^ 
are conserved at every node. It is also worth noting that larger values of C, and 
therefore a better efficiency of the the mean-field approximation, are found for large 
temperatures. For instance, at T = 5 we get C = 0.63; as the temperature is lowered 
the correlation decreases displaying a possible minimum around Tc. Critical effects 
apart, the mean-field approach seems to provide good estimates especially for large 
temperatures: indeed for T approximately larger than 3.5 one has C > 0.5 and we 
checked that at T = 500, C is close to 0.94. We also underline that C = 1 means that 
local conductivities are purely local quantities, independent of the neighborhood. 

We remark that the links crossed by large currents may be sharply identified 
within this approximation, which, for example, captures the link evidenced in Fig. [6] 
although its intermediate coupling value, it carries a large current. In other terms, 
our approximation is able to locate the regions characterized by high currents which, 
in realistic realizations, could lead to a failure of the link itself. 



6. Conclusions and Perspectives 

In this work we addressed the general problem of a spin model on arbitrary discrete 
structures with a microcanonical dynamics, focusing on the conduction properties at 
microscopic scales. First, we introduced a consistent definition of local currents. Then 
we applied such general definition in the case of disordered ferromagnetic Ising model 
on a cylindrical structure, where boundaries are coupled with thermostats at different 
fixed temperatures. We highlighted that the local microscopic currents depend non 
trivially on the whole distribution of quenched couplings and a consistent definition 
of local conductivity has been introduced, at least in regime of linear response, i.e. 
when the temperature gradient AT is small. 
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In spite of the aforementioned dependence of local quantities (currents, 
temperatures and conductivities) on the whole disorder arrangement, numerical results 
evidenced a simple correlation between local conductivities and the pertaining local 
couplings. These correlations suggested the development of an approximated mean- 
field like approach, where local conductivities depend on local couplings only, and 
the conduction properties (i.e. the currents) can be easily evaluated by solving 
the discretized linear Fourier equations with suitable border conditions. Such 
approximation is especially effective at large temperatures. With respect to previous 
work on the cylindrical spin model, the strict requirement on the smallness of AT is 
due precisely to the fact that in the present case we focus on the microscopic aspects. 

This pattern of results naturally indicates preferential lines of future 
developments. First of all, the definition of currents holds on a general geometrical 
substrate. In a complex topological structure, some aspects of the framework we have 
worked out are expected to be still valid, in particular the presence of a linear response 
regime and the definition of local conductivities whose values may depend, however, 
in a non trivial way on the topology of the underlying structure. 

Moreover, an important point to investigate concerns the extent of correlation 
length, as a function of temperature, degree of disorder and topology. In other terms, 
it would be interesting to determine the length of the radius such that the external 
pattern constitutes a practically uniform background, without any influence on the 
local currents and conductivities. Also the response of local conductivity to small 
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coupling perturbations would be in order. This kind of problems, at finite and low 

temperatures, has obviously to do with the correct definition of the mesoscopic scale 
we spoke about. At T :» 1, analytical and numerical estimates say that this radius 
possibly reduces to one link, which means that the conductivity is a local property. 

We also infer the interest of a more precise characterization of the behaviour of 
the system both at low temperatures and for T « Tc, where the linear approximation 
has some difficulties, possibly not only of numerical nature. Finally, the whole set of 
results clearly indicates that we have an extremely flexible calculation engine ready 
for more complex realizations of disorder, including the case of null and negative 
couplings, i.e. diluted and spin-glass systems. 
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